Intersections of Adversity & Neurodiversity: Adverse Childhood Experiences’ Association to Mental Health & the Buffering Role of Flourishing
This study extends analyses conducted by Bethell et al. (2023), who explored the role of family resilience and connection in promoting flourishing among U.S. children facing adversity. Our study applies a similar framework but focuses on neurodivergent children, examining how adverse childhood experiences (ACEs) influence mental health outcomes and how flourishing behaviors—such as curiosity, emotional control, and task persistence—may buffer these effects. This study was conducted by the Slopen Laboratory at the Harvard T.H. Chan School of Public Health.
Author: Adrián Medina
Date: October 5, 2024
Project Overview
This repository contains files relevant to a study investigating the relationship between adverse childhood experiences (ACEs) and mental health outcomes in children with neurodevelopmental disorders (NDDs), such as Autism Spectrum Disorder (ASD) and ADHD. The study specifically examines the moderating role of child flourishing behaviors—curiosity in learning, emotional control, and task persistence—on the effects of ACEs on mental health outcomes.
Study Design:
The study utilizes data from the National Survey of Children’s Health (NSCH), focusing on children aged 6-17 with reported neurodevelopmental issues (N = 44,776, MAge = 12.2). This cross-sectional analysis examines ACE exposure and its impact on three primary mental health outcomes: anxiety, depression, and behavioral issues. The moderating role of child flourishing behaviors is assessed using interaction terms in logistic regression models.
Data Collection:
Key variables collected include:
Adverse Childhood Experiences (ACEs): Derived from parent-reported experiences such as exposure to violence, household dysfunction, and discrimination.
Mental Health Outcomes: Parent-reported diagnoses of anxiety, depression, and behavioral problems, confirmed by health care providers or educators.
Child Flourishing Metrics: Based on parent-reported behaviors such as curiosity in learning, emotional control, and task persistence, which were used to develop a Child Flourishing Index (CFI).
Analytic Approach:
The primary analytic approach involved logistic regression models to assess the association between ACEs and mental health outcomes, with interaction terms to explore the moderating effect of child flourishing. Key covariates include age, sex, race/ethnicity, and socioeconomic status.
Goals:
The study aims to:
Examine the dose-response relationship between the number of ACEs and risks of mental health challenges in neurodivergent children.
Identify the protective role of child flourishing behaviors, offering insights into resilience-promoting interventions for this vulnerable population.
Tip
For detailed information on the dataset variables, refer to the NSCH Data Dictionary.
Set up the R environment by configuring the CRAN repository, installing essential packages, setting the base path, and loading the primary dataset into the local environment.
Expand Code
# Set CRAN repository for consistent package installationoptions(repos =c(CRAN ="http://cran.rstudio.com/"))# Install and load the pacman package for efficient package managementif (!require(pacman)) install.packages("pacman")library(pacman)# Use p_load function to install (if necessary) and load packagesp_load(dplyr, tidyr, tidyverse, knitr, ggplot2, Hmisc, broom, stats, kableExtra, webshot2, sjPlot)# Specify the 'base_path' where you can find your data files, ASSUMING they're in the same directory, & set it as WDbase_path <-"~/GitHub/Adversity-Neurodiversity/data_files"setwd(base_path)# Load primary data fileneurodivergent_data <-read_csv("neurodivergent_data.csv")
1
See details under the important callout below!
Data Extraction & Filtering (Archived Code)
The file being used for these analyses is a subset of a “master” file containing the compiled contents of the data releases from the National Survey on Child Health’s 2016-2022 cycles. The master file is nearly 500 MB (0.5 GB) in size, making it cumbersome to process in real-time. For efficiency and ease of analysis, I am using a smaller subset that contains only the relevant data. This subset was generated using the data extraction and filtering process detailed in the archived code below, in an effort to maintain transparency and ensure reproducibility of the workflow.
Expand Code
# Define the path to the datasetdata_path <-"~/Downloads/Data2_2016to2022.csv"Data2_2016to2022 <-read.csv(data_path)# Prepare project data: Select specific variables, create 'neurodivergent' variable, & filter data based on age and neurodivergent statusneurodivergent_data <- Data2_2016to2022 %>%select(year, fpl, fpl_mean, sc_age_years, sc_hispanic_r, sc_race_r, sc_sex, higrade_tvis, k2q35a, k2q35c, k2q38a, k2q38c, k2q30a, k2q30c, k2q36a, k2q36c, k2q60a, k2q60c, k2q37a, k2q37c, downsyn, downsyn_desc, k2q31a, k2q31c, k2q61a, cerpals_desc, k2q33a, k2q33b, k2q33c, k2q32a, k2q32b, k2q32c, k2q34a, k2q34b, k2q34c, ace1, ace3, ace4, ace5, ace6, ace7, ace8, ace9, ace10, ace11, ace12, k6q71_r, k7q85_r, k7q84_r) %>%# Creating 'neurodivergent' variable as NSCH does not explicitly inquire about neurodivergent status.# Using a list of diagnoses provided by NSCH to define neurodivergent individuals.mutate(neurodivergent =if_else(k2q35a ==1| k2q38a ==1| k2q36a ==1| k2q60a ==1| k2q37a ==1| k2q30a ==1| downsyn ==1| k2q31a ==1| k2q61a ==1, 1, 0)) %>%# Filter data to include only individuals aged 6 or above and identified as neurodivergentfilter(sc_age_years >=6, neurodivergent ==1)# Resulting filtered data to be used for further analysiswrite.csv(neurodivergent_data, "neurodivergent_data.csv")
Analytic Data Preparation
Set up recoding of variables by transforming ‘ace1’ into a dichotomous variable, calculating the Childhood Flourishing Index (CFI), and categorizing total ACE counts. Recode education, income, sex, race/ethnicity, and CFI into categorical factors for analysis.
Expand Code
# Recode 'ace1' to dichotomous variable based on analytic specificationsneurodivergent_data$ace1_recode <-ifelse(neurodivergent_data$ace1 %in%c(2, 3, 4), 1, ifelse(neurodivergent_data$ace1 ==1, 0, NA))# Recode ACE counts and calculate Childhood Flourishing Index (CFI)neurodivergent_data <- neurodivergent_data %>%mutate(# Total ACE count for each child excluding missing valuesace_total =rowSums(select(., ace1_recode, ace3:ace12) ==1, na.rm =TRUE),# Categorize total ACEs for further analysisace_counts =cut(ace_total, breaks =c(-1, 0, 1, 3, Inf), labels =c('0', '1', '2-3', '4+'), right =TRUE),ace_counts =factor(ace_counts, levels =c("0", "1", "2-3", "4+")), # Recreate 'ace_counts' factor# Calculate the Childhood Flourishing Index (CFI)CFI =rowSums(select(., k6q71_r, k7q85_r, k7q84_r) ==1, na.rm =TRUE),# Recode CFI into dichotomous variableCFI_dich =case_when( CFI ==3~"Flourishing", CFI %in%0:2~"Not Flourishing",TRUE~NA_character_ ),CFI_dich =factor(CFI_dich, levels =c("Not Flourishing", "Flourishing")) # Recreate 'CFI_dich' factor )# Recode binary outcomes for anxiety, depression, and behavioral issues into factorsneurodivergent_data <- neurodivergent_data %>%mutate(# Recode Anxiety (k2q33b)Anxiety =factor(k2q33b, levels =c(2, 1), labels =c("No", "Yes")),# Recode Depression (k2q32b)Depression =factor(k2q32b, levels =c(2, 1), labels =c("No", "Yes")),# Recode Behavioral Issues (k2q34b)Behavioral_Issues =factor(k2q34b, levels =c(2, 1), labels =c("No", "Yes")) )# Recode education and income categoriesneurodivergent_data <- neurodivergent_data %>%mutate(# Recode education into categorical factorshighgrade_tvis_cat =case_when( higrade_tvis ==1~"Less than high school", higrade_tvis ==2~"High school (including vocational)", higrade_tvis ==3~"Some college or associate degree", higrade_tvis ==4~"College degree or higher",TRUE~NA_character_ ),highgrade_tvis_cat =factor(highgrade_tvis_cat, levels =c("Less than high school", "High school (including vocational)", "Some college or associate degree", "College degree or higher")),# Recode federal poverty level categories into factorsfpl_cat =case_when( fpl %in%c(50:99) | fpl_mean <100~"Less than 100%", fpl %in%c(100:199) | fpl_mean <200~"100% to 199%", fpl %in%c(200:299) | fpl_mean <300~"200% to 299%", fpl %in%c(300:399) | fpl_mean <400~"300% to 399%", fpl %in%c(400:999) | fpl_mean <999~"400% or greater",TRUE~NA_character_ ),fpl_cat =factor(fpl_cat, levels =c("Less than 100%", "100% to 199%", "200% to 299%", "300% to 399%", "400% or greater")) )# Recode sex and race categoriesneurodivergent_data <- neurodivergent_data %>%mutate(# Recode sex into categorical factorssc_sex_cat =case_when( sc_sex ==1~"Male", sc_sex ==2~"Female",TRUE~NA_character_ ),sc_sex_cat =factor(sc_sex_cat, levels =c("Male", "Female")),# Recode race and ethnicityrace =case_when( sc_race_r ==1~"White", sc_race_r ==2~"Black", sc_race_r ==3~"American Indian or Alaska Native", sc_race_r %in%4:5~"Asian, Native Hawaiian, or Pacific Islander", sc_race_r %in%6:7~"Other or mixed race",TRUE~NA_character_ ),ethnicity =case_when( sc_hispanic_r ==1~"Hispanic/Latino", sc_hispanic_r ==2~"Non-Hispanic/Latino",TRUE~NA_character_ ),# Recreate 'sc_race_cat' from race and ethnicity combinationssc_race_cat =case_when( race %in%c("White") & ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic White", race %in%c("Black") & ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Black or African American", race %in%c("American Indian or Alaska Native") & ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic American Indian or Alaska Native", race %in%c("Asian, Native Hawaiian or Pacific Islander") & ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Asian, Native Hawaiian or Pacific Islander", race %in%c("Other or mixed race") & ethnicity =="Non-Hispanic/Latino"~"Non-Hispanic Other or mixed race", ethnicity =="Hispanic/Latino"~"Hispanic or Latino of any race",TRUE~NA_character_ ),sc_race_cat =factor(sc_race_cat, levels =c("Non-Hispanic White", "Non-Hispanic Black or African American", "Non-Hispanic American Indian or Alaska Native", "Non-Hispanic Asian, Native Hawaiian or Pacific Islander", "Non-Hispanic Other or mixed race", "Hispanic or Latino of any race")) )
1
ace1 is originally coded by NSCH as a frequency measure, but it’s needed as indicative variable of presence like the other ace# variables.
2
Part of this includes the calculating of the Childhood Flourishing Index (CFI) based on Bethell et al. (2019) criteria
Frequencies & Descriptive Statistics
Expand Code
# Creating frequency tables for specific diagnostic outcomes to assess their distribution in the datasetfreq_k2q33b <-table(neurodivergent_data$k2q33b ==1, useNA="ifany")freq_k2q32b <-table(neurodivergent_data$k2q32b ==1, useNA="ifany")freq_k2q34b <-table(neurodivergent_data$k2q34b ==1, useNA="ifany")# Combine frequency tables into a single data frame for easy viewingfreq_table <-data.frame(Diagnosis =c("Anxiety", "Depression", "Behavioral Problems"),Count =c(freq_k2q33b["TRUE"], freq_k2q32b["TRUE"], freq_k2q34b["TRUE"]))# Use knitr::kable to format the frequency table for better readability in the output documentknitr::kable(freq_table, caption ="Frequency of Specific Diagnostic Outcomes")
Frequency of Specific Diagnostic Outcomes
Diagnosis
Count
Anxiety
13725
Depression
6289
Behavioral Problems
13535
Expand Code
# Count the number of subjects with a CFI score of 3 using a logical condition, as this is considered 'Flourishing.'num_subjects_cfi_3 <-sum(neurodivergent_data$CFI ==3, na.rm =TRUE)cat("Number of Subjects Flourishing (CFI Score of 3):", num_subjects_cfi_3, "\n")
Number of Subjects Flourishing (CFI Score of 3): 3265
Expand Code
# Histogram of ages for neurodivergent participantsggplot(neurodivergent_data, aes(x = sc_age_years)) +geom_histogram(binwidth =1, fill ="blue", color ="black") +labs(title ="Age Distribution of Neurodivergent Participants", x ="Age (Years)", y ="Count") +theme_minimal()
Expand Code
# Bar plot of ACE counts categoriesggplot(neurodivergent_data, aes(x =factor(ace_counts), fill =factor(ace_counts))) +geom_bar() +geom_text(stat ='count', aes(label = ..count..), position =position_stack(vjust =0.5), color ="white") +labs(x ="ACE Counts Category", y ="Number of Children", fill ="ACE Category",title ="Distribution of Children by ACE Counts") +scale_x_discrete(labels =c('1'="0 ACEs", '2'="1 ACE", '3'="2-3 ACEs", '4'="4+ ACEs")) +theme_minimal() +theme(axis.text.x =element_text(angle =45, hjust =1))
Expand Code
# Bar plot of CFI counts categoriesggplot(neurodivergent_data, aes(x =factor(CFI), fill =factor(CFI))) +geom_bar() +geom_text(stat ='count', aes(label = ..count..), vjust =-0.5, color ="black") +labs(x ="Childhood Flourishing Index (CFI) Score",y ="Number of Children",fill ="CFI Score",title ="Distribution of Children by CFI Score") +scale_x_discrete(labels =c('0'="0", '1'="1", '2'="2", '3'="3")) +theme_minimal() +theme(axis.text.x =element_text(angle =45, hjust =1))
Expand Code
# Bar plot of CFI_dich counts categoriesggplot(neurodivergent_data, aes(x = CFI_dich)) +geom_bar(aes(fill = CFI_dich)) +geom_text(stat ='count', aes(label = ..count..), vjust =-0.5) +scale_fill_manual(values =c("0"="blue", "1"="green")) +labs(title ="Distribution of CFI Dichotomous Variable",x ="Child Flourishing Index (Dichotomous)",y ="Count") +theme_minimal() +theme(legend.title =element_blank()) # Remove the legend title if not needed
The tab_model function from the sjPlot package automatically exponentiates the logistic regression coefficients into Odds Ratios when displaying models with a binomial family (logistic regression).
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Expand Code
# Plot for the Anxiety model with x-axis limits between 0.1 and 5plot_model(model_anx, type ="est", show.values =TRUE, show.p =TRUE, ci.lvl =0.95, title ="Odds Ratios - Anxiety", axis.lim =c(0.1, 5))
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Expand Code
# Plot for the Depression model with x-axis limits between 0.1 and 5plot_model(model_dep, type ="est", show.values =TRUE, show.p =TRUE, ci.lvl =0.95, title ="Odds Ratios - Depression", axis.lim =c(0.1, 5))# Plot for the Behavioral Issues model with x-axis limits between 0.1 and 5plot_model(model_beh, type ="est", show.values =TRUE, show.p =TRUE, ci.lvl =0.95, title ="Odds Ratios - Behavioral Issues", axis.lim =c(0.1, 5))
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Anxiety
Depression
Behavioral Issues
Predictors
Odds Ratios
CI
p
Odds Ratios
CI
p
Odds Ratios
CI
p
(Intercept)
11.02
6.82 – 18.33
<0.001
6.93
4.08 – 11.99
<0.001
49.42
34.07 – 72.76
<0.001
ace counts [1]
1.29
1.09 – 1.53
0.003
1.24
1.02 – 1.51
0.033
1.14
0.99 – 1.31
0.063
ace_counts2-3
1.52
1.28 – 1.81
<0.001
1.47
1.22 – 1.78
<0.001
1.12
0.98 – 1.29
0.095
ace counts [4+]
1.88
1.54 – 2.29
<0.001
1.80
1.47 – 2.20
<0.001
1.36
1.17 – 1.59
<0.001
sc race cat [Non-Hispanic
Black or African
American]
0.76
0.59 – 1.01
0.053
0.94
0.72 – 1.25
0.658
1.03
0.86 – 1.24
0.742
sc race cat [Non-Hispanic
American Indian or Alaska
Native]
0.52
0.31 – 0.93
0.018
1.35
0.68 – 3.06
0.430
0.91
0.57 – 1.52
0.707
sc race cat [Non-Hispanic
Other or mixed race]
0.89
0.71 – 1.14
0.349
1.01
0.80 – 1.29
0.917
0.93
0.78 – 1.10
0.387
sc race cat [Hispanic or
Latino of any race]
0.70
0.59 – 0.84
<0.001
0.84
0.69 – 1.02
0.071
0.90
0.78 – 1.05
0.169
sc sex cat [Female]
1.35
1.19 – 1.53
<0.001
1.42
1.25 – 1.61
<0.001
1.21
1.09 – 1.35
0.001
fpl cat [100% to 199%]
0.95
0.76 – 1.19
0.670
0.88
0.70 – 1.10
0.273
0.91
0.77 – 1.09
0.310
fpl cat [200% to 299%]
1.03
0.81 – 1.30
0.799
0.82
0.65 – 1.04
0.109
0.83
0.69 – 0.99
0.039
fpl cat [300% to 399%]
0.84
0.66 – 1.07
0.152
0.78
0.61 – 1.01
0.057
0.69
0.57 – 0.83
<0.001
fpl cat [400% or greater]
0.91
0.72 – 1.15
0.443
0.68
0.54 – 0.86
0.002
0.71
0.59 – 0.84
<0.001
highgrade tvis cat [High
school (including
vocational)]
1.12
0.72 – 1.69
0.605
1.10
0.72 – 1.65
0.651
1.00
0.72 – 1.36
0.981
highgrade tvis cat [Some
college or associate
degree]
1.11
0.72 – 1.65
0.627
1.14
0.74 – 1.68
0.539
0.98
0.71 – 1.33
0.887
highgrade tvis cat
[College degree or
higher]
1.22
0.79 – 1.83
0.353
1.04
0.68 – 1.54
0.864
0.86
0.62 – 1.16
0.330
sc age years
0.97
0.95 – 0.99
0.001
0.96
0.93 – 0.98
<0.001
0.85
0.84 – 0.87
<0.001
Observations
14612
7372
15339
R2 Tjur
0.007
0.018
0.040
Expand Code
# Logistic regression model with CFI_dich as a moderator for Anxietymodel_anx_mod <-glm(Anxiety ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_cat + highgrade_tvis_cat + sc_age_years,data = neurodivergent_data, family = binomial)# Logistic regression model with CFI_dich as a moderator for Depressionmodel_dep_mod <-glm(Depression ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_cat + highgrade_tvis_cat + sc_age_years,data = neurodivergent_data, family = binomial)# Logistic regression model with CFI_dich as a moderator for Behavioral Issuesmodel_beh_mod <-glm(Behavioral_Issues ~ ace_counts * CFI_dich + sc_race_cat + sc_sex_cat + fpl_cat + highgrade_tvis_cat + sc_age_years,data = neurodivergent_data, family = binomial)# Display the models with tab_modeltab_model(model_anx_mod, model_dep_mod, model_beh_mod, show.ci =0.95, transform =NULL)
1
Since we included interaction terms, we want to see the raw coefficients instead. You can set transform = NULL in the tab_model function to ensure the values aren’t exponetiated into Odds Ratios.
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Profiled confidence intervals may take longer time to compute.
Use `ci_method="wald"` for faster computation of CIs.
Expand Code
# Plot marginal effects for the Anxiety modelplot_model(model_anx_mod, type ="int", title ="Interaction: ace_counts * CFI_dich (Anxiety)")# Plot marginal effects for the Depression modelplot_model(model_dep_mod, type ="int", title ="Interaction: ace_counts * CFI_dich (Depression)")# Plot marginal effects for the Behavioral Issues modelplot_model(model_beh_mod, type ="int", title ="Interaction: ace_counts * CFI_dich (Behavioral Issues)")
Anxiety
Depression
Behavioral Issues
Predictors
Log-Odds
CI
p
Log-Odds
CI
p
Log-Odds
CI
p
(Intercept)
2.45
1.97 – 2.96
<0.001
1.96
1.43 – 2.51
<0.001
3.95
3.58 – 4.34
<0.001
ace counts [1]
0.26
0.08 – 0.43
0.004
0.20
-0.00 – 0.40
0.053
0.14
-0.00 – 0.28
0.053
ace_counts2-3
0.41
0.23 – 0.58
<0.001
0.34
0.15 – 0.53
0.001
0.12
-0.02 – 0.26
0.090
ace counts [4+]
0.59
0.39 – 0.79
<0.001
0.54
0.34 – 0.75
<0.001
0.31
0.15 – 0.46
<0.001
CFI dich [Flourishing]
-0.96
-1.43 – -0.46
<0.001
-1.43
-2.23 – -0.66
<0.001
-0.76
-1.31 – -0.17
0.009
sc race cat [Non-Hispanic
Black or African
American]
-0.28
-0.54 – 0.00
0.047
-0.06
-0.33 – 0.22
0.674
0.03
-0.15 – 0.22
0.735
sc race cat [Non-Hispanic
American Indian or Alaska
Native]
-0.66
-1.18 – -0.07
0.018
0.28
-0.40 – 1.11
0.452
-0.11
-0.57 – 0.41
0.671
sc race cat [Non-Hispanic
Other or mixed race]
-0.13
-0.36 – 0.11
0.287
0.01
-0.23 – 0.25
0.953
-0.07
-0.24 – 0.10
0.412
sc race cat [Hispanic or
Latino of any race]
-0.35
-0.53 – -0.16
<0.001
-0.18
-0.37 – 0.02
0.076
-0.10
-0.25 – 0.04
0.165
sc sex cat [Female]
0.31
0.18 – 0.43
<0.001
0.35
0.23 – 0.48
<0.001
0.19
0.08 – 0.30
0.001
fpl cat [100% to 199%]
-0.05
-0.27 – 0.18
0.692
-0.12
-0.35 – 0.10
0.296
-0.11
-0.29 – 0.06
0.211
fpl cat [200% to 299%]
0.03
-0.21 – 0.27
0.792
-0.19
-0.43 – 0.04
0.112
-0.21
-0.39 – -0.03
0.024
fpl cat [300% to 399%]
-0.17
-0.41 – 0.07
0.174
-0.25
-0.51 – 0.00
0.055
-0.40
-0.59 – -0.21
<0.001
fpl cat [400% or greater]
-0.08
-0.32 – 0.15
0.487
-0.38
-0.61 – -0.14
0.002
-0.37
-0.55 – -0.19
<0.001
highgrade tvis cat [High
school (including
vocational)]
0.08
-0.36 – 0.50
0.700
0.10
-0.33 – 0.50
0.650
-0.02
-0.36 – 0.29
0.880
highgrade tvis cat [Some
college or associate
degree]